package scu.maqiang.homogenization;

import scu.maqiang.fes.BVPType;
import scu.maqiang.fes.FES2T32;
import scu.maqiang.fes.Materials;
import scu.maqiang.mesh.MatrixFunc;
import scu.maqiang.mesh.Mesh2T3;
import scu.maqiang.mesh.ScalarFunc;
import scu.maqiang.numeric.MVO;

public class ElasticityEigenValue2DT3 {

	public static void main(String[] args) {
		long begin = System.nanoTime();
		int N = 8;
		double ep = 1.0 / N;
		
		Mesh2T3 cellMesh = new Mesh2T3().square(32, 32);
		ScalarFunc regionCellFunc = (xy, llable, param) -> {
            return (Math.abs(xy[0] - 0.5) < 0.25 && Math.abs(xy[1] - 0.5) < 0.25)? 1.0:-1.0;
        };
	    cellMesh.setDomainLabel(regionCellFunc, null, 10);
	    FES2T32 fsCell = new FES2T32(cellMesh);
	    double[][] FOCS = new double[3][fsCell.GetNdof()];
		double[][] SOCS = new double[6][fsCell.GetNdof()];
		double[][] D = new double[3][3];
		SOTS2DT362 ss = new SOTS2DT362(fsCell, null);
		
		double E1 = 1.5e4, E2 = 2.01e6;
		double Nu1 = 0.3, Nu2 = 0.25;
		double[][] D1 = Materials.ElasticConstituteMatrix(E1, Nu1, BVPType.PLANESTRAIN);
		double[][] D2 = Materials.ElasticConstituteMatrix(E2, Nu2, BVPType.PLANESTRAIN);
		MatrixFunc mf = (xyz, label, param) -> label == 10? D1 : D2;
		ScalarFunc rhoF = (xyz, label, param) -> label == 10? 0.5 : 1.0;
		double[] homoRho = {0.0};
		
		ss.elasticityEigenValueCellSolution(mf, new ScalarFunc[] {rhoF}, null, BVPType.CONSTITUTE_MATRIX_COMMON, FOCS, D, homoRho, SOCS);
		
		System.out.println("min: " + MVO.min(FOCS[0]) + "\tmax: " + MVO.max(FOCS[0]));
		System.out.println("min: " + MVO.min(FOCS[1]) + "\tmax: " + MVO.max(FOCS[1]));
		System.out.println("min: " + MVO.min(FOCS[2]) + "\tmax: " + MVO.max(FOCS[2]));
		System.out.println("Homogenized Elasticity: ");
		System.out.println(MVO.toString(D));
		System.out.println("Homogenized Density:");
		System.out.println(homoRho[0]);
		System.out.println("min: " + MVO.min(SOCS[0]) + "\tmax: " + MVO.max(SOCS[0]));
		System.out.println("min: " + MVO.min(SOCS[1]) + "\tmax: " + MVO.max(SOCS[1]));
		System.out.println("min: " + MVO.min(SOCS[2]) + "\tmax: " + MVO.max(SOCS[2]));
		System.out.println("min: " + MVO.min(SOCS[3]) + "\tmax: " + MVO.max(SOCS[3]));
		System.out.println("min: " + MVO.min(SOCS[4]) + "\tmax: " + MVO.max(SOCS[4]));
		System.out.println("min: " + MVO.min(SOCS[5]) + "\tmax: " + MVO.max(SOCS[5]));
	}

}
